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ABSTRACT 

The development of a numerical 1000=-mb prognosis by means 
of a 500 to 1000—mb thickness forecast is investigated. The 
thickness forecast is based, in part,on the 500=-mb barotropic 
prognosis as developed by Fleet Numerical Weather Facility, 
Monterey. 

Forecasts were prepared for three days during May 1966, 
with various amounts of smoothing a the component fields. 
The model showed promise of producing a rapid, usable prognosis; 
however further testing is needed with the probable addition of 
some climatic constraints. It is also hoped that the yet untried 
diabatic heating and surface friction terms will provide some 


further improvement. 


Section 


DUDLEY KNOX LIBRARY 
NAVAL POSTGRADUATE SCHOOL 
MONTEREY, CA 93943-5101 


TABLE OF COM PENTsS 


Iptroduetion 

Development of the Prognostic Model 
Development of the Heating Term 

Development of the Terrain and Friction Terms 


Programming the Prognostic Equations for the 
Control Data 1604 Computer 


Results and Conclusions 
Acknowledgements 


Bibliography 


Page 
gs 
12 
16 


ad 


Ze 
Zc 
29 


30 





10. 


ii 


are 


13. 


1h. 


aS 


LIST OF ILLUSTRATIONS 


Two Layer Division of the Atmosphere 
between 1000—mb and 500-—mb 


Plot of Bowen Ratio vs Sea-Air Temperature 
Difference 


500-1000 mb Thickness Analysis 0024 
20 May 1966 


24-Hour Thickness Prognosis from 002 
20 May 1966 


2ehe-Hour Thickness Prognosis from 002 
20 May 1966 (Decreased Smoothing) 


500-1000 mb Thickness Analysis 002 
21 May 1966 


1000=mb Height Analysis 002 
20 May 1966 


2y-Hour 1000—mb Height Frognosis from 002% 
20 May 1966 


500-1000 mb Thickness Analysis 122 
23 May 1966 


2ey—-Hour Thickness Prognosis from 1224 
23 May 1966 


500-1000 mb Thickness Analysis 122 
2, May 1966 


1000—mb Height Analysis 122 
23 May 1966 


LheHour 1000—-mb Height Prognosis from 1224 
23 May 1966 


500-1000 mb Thickness Analysis 122 
29 May 1966 


24—-Hour Thickness Prognosis from 122 
29 May 1966 


Page 


So. 


Se 


33 


34 


5): 


36 


Sy 


BS 


39 


1.0 


41 


he 


43 


hh 


45 


16. 


17. 


18. 


19. 


20. 


‘el ae 


ENWE 24—Hour Thickness Prognosis from 
122 29 May 1966 


500-1000 mb Thickness Analysis 122 
30 May 1966 


1000-mb Height Analysis 122 29 May 1966 


2Ly-Hour 1000—-mb Height Prognosis from 122 
29 May 1966 


ENWE 2heHour 1000-mb Height Analysis from 
122 29 May 1966 


1000-mb Height Analysis 122 30 May 1966 


46 


47 
48 


49 


50 
Sn 


8 


< 


RMSE 


TABLE OF SYMBOLS AND ABBREVIATIONS 


Specific heat of dry air at constant pressure 
Specific heat of dry air at constant volume 
Two dimensional del operator on constant pressure surface 
Coriolis parameter 

Fleet Numerical Weather Facility, Monterey, Califomia 
Acceleration of gravity 

Thickness of layer 1000 mb to 500 mb 
Thickness of layer 750 to 500 mb 

Thickness of layer 1000 mb to 500 mb 
Jacobian 

Latent heat of evaporation 

millibar 

Absolute voriticity 

Pressure 

Latitude (degrees) 

Change of thickness due to diabatic heating 
Heating rate 

Flux of heat through air/ocean interface 
Gas constant per gram of dry air 
Root-mean~square error 

Density 


Time 


Temperature (°K) 

Temperature of air at the surface 
Vector wind 

Vertical velocity 

Terrain induced vertical velocity 
Height of constant pressure surface 


Relative voriticity 


gs 


gsl0 


10 


SUBSCRIPTS 
air at surface 
previous estimate 
geostrophic 
1000-mb geostrophic 
sea surface 
terrain 
900-mb 
625=mb 
750-mb 
875-mb 


1000-mb 








leo Tiittreauction 

With the development of the high speed digital computer, the 
problem of weather prognosis has evolved from one of qualitative 
reasoning by an experienced forecaster to one of objective 
computations based on mathematical equations representing the 
motions of the atmosphere. At this time, daily numerical analysis 
and prognostic charts are being produced at Fleet Numerical Weather 
Facility, Monterey (FNWF), and rapidly disseminated to fleet units 
around the world. The entire field of meteorology is in a rapid state 
of change with more and more products being produced numerically. 

The most successful numerical prognosis has been at the 500 mb 
level. Ifa successful numerical isobaric thickness forecast could 
be developed, the success of the 500emb forecast could then be 
reflected to higher and lower levels. The objective of this research 
is to obtain a dependable 24-hour numerical forecast of thickness 
for the layer between 500-mb and 1000-mb, and then using FNWF's 
500-mb barotropic forecast obtain a 1000-mb height prognosis. 

The thickness of an air layer at any given time and place is 
dependent upon the mean temperature of the air column. If we look 
at thickness change as a result of the advection of thickness lines 
only, it appears that there is no’ heat loss or gain as the air moves 
across the surface of the earth nor are there vertical motions. This 


is not true, of course, Since the air gains or loses heat through 


zal 


radiative, dynamical, and conductive processes, and through latent 
heat exchange resulting from phase changes of the water constituent. 

In the model used in this study, in addition to advection, an 
empirical heating term is developed to acéount for heating or cooling 
by evaporation or condensation and for turbulent transfer of 
sensible heat over ocean surfaces. 

To account for the dynamic processes, vertical velocity is 
assumed to be a second order function of pressure. The 1000-mb 
vertical velocity is composed of a terrain upslope term and a 
friction term. Heat gained or lost by radiative processes was 
ignored as neglibible — the 24-hour period. 

2. Development of the Prognostic Model 

A two layer model is used in the development of the prognostic 
equations. The lower troposphere is divided into two layers as. 
illustrated in Fig. l. 

By use of the equation of state and the hydrostatic equation, 
the thermodynamic equation may be reduced to the following form 


&(%) +V-VEE +0 (BF +22) = - 80, 4) 


i 


Replacing the vertical derivatives above by finite difference and 


Simplifying leads to 


ab Yh + 0 (Sph—bzhl)- BBQ 


oP 


IL 


Where AZ=h , the thickness. 
In the layer 1000 to 500 mb it has been noted that @ can 


generally be approximated by a second order equation in pressure. 


W= A(xyt) (tat) + DOKI ye) + We (3) 


By means of equation(2)it may be shown that 


VW = Yea + yb + We (4) 
We = Vy ® + Vy. © + Wro (5) 
OW) - _AQ — De ol) = =. (6) 
(3), Pe Ps j (SP) Ps 


Taking the vorticity equation at 500 mb in the form 


a V4, = W, (3), 


and substituting from equation ‘(6) gives 


= -20 J ( Qt +\. -\/™.) (7) 


Equation (2) is now evaluated for layer 2 giving 


a VW -Vh, + MB hy — Vhs rh Bos core 


The above equation may -be:written 


= £4 -B-YT,) 
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Where 


0, =(&, h, — (b ba b)) 


Combing equations (5) and (8) and simplifying leads to 


ele aN. 3 uy), \-& B (36 )i9) 
7 La bGeh Co“ Se WVha)- We) -$ f(g? NV, 
With the geostrophic assumption for relative vorticity 
Y= H7* oz 
at EV ot 
it follows that 


4h Ve = =, vs, ai ™ (SP 2h 
LV" Qz = — V, Viele + (35), (ata) 


Subtracting equation (11) from (10 yields 


Ye 7 “dhs =," ‘V4,.-V, V4, +34), amg ge) (12) 


dt 


However, 


Qu) = (We = We) » (Ag- We) = -L (asd) 


1 ae 


ee aaa 


thus equation (12 becomes 


M7 Qh =\V,V4%,-V V4, +e 8 ot (ath) (13) 
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Substituting equations (7) and (9) into equation (13) and expressing 


advective terms in Jacobian form we obtain the following: 
at TEAS a om J(D. Yo) J (0, 1) "Tp. GJ(0,, D,) 


bt GO: DM) + FV 2? — 165 G RAP ¢ 
TERS TE Mee HF I Date) + F VG — HF CRO 


fleet. 


(14) 


where 


WV ha = Be J (Z 10 29) = 3] J (Do,D,) W, V0 J(Do,.) 


and 
F= (Mot34,) » G= (4 - ae) 
74, Oa 

The prognostic thickness equation is thus determined fa the lower 
layer. This equation is of the Helmholtz type, the solution of 
which may readily be obtained by use of a digital computer. 

We will riow obtain a prognostic equation for the thickness 
of the upper layer. 

Evaluating equation (2) for the upper level centered at 625 mb 


results in the following: 
Yh, — . . - A @ 
at Y. Vb, ‘ GW, aCe Q, (AES) 
where os C ( i h 
— a) ae — ee 2. ) 
( Cet h, as r 
Substituting from equations (4), (7) and (9) and simplifying, we 


obtain 


OW), = “See ham x St Aa (SE) REM) (16) 


[S_ fA ° 
Ls Aas = O41 
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Substitution of equation (16) into equation (15), and expressing 
in Jacobian form provides the thickness tendency (prognostic) 


equation for the upper layer of the model. 


Hs (D;,0;) )+ Xs 2 FD» D,) D, He HS Sh Seite VOR (17) 


* Taig Pat J Osi, ere 26,48 2 GWao 


where 
VVh, = 35 J (Orde) 


Ill. Development of the heating term. 
Initially it will be assumed that the heating rate is expressed 


by an equation of the following format 


Q “AS MMT.-T)(4£) (18) 


Over land the heating term is quite small, and we shall 
consider it as negligible. Over water, excluding radiation, the 
main contributions to the heating term are due to sensible and 
latent heat exchange between the sea and the atmosphere. Itis 
of course realized that latent heat gained through evaporation may 
not be utilized by the atmosphere at the exact instant and position 
that it is obtained. However, in the relationship derived, an 
empirical adjustment constant is available should the values 


calculated be determined to be somewhat too small or large. 
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According to Laevastu [3] the heat flux from the ocean to 
the atmosphere due to sensible heat transfer and evaporation or 


condensation may be expressed as follows: 


Qea = 34 (0.26 +0.71Ve)(Ts —TVa) (aace. (19) 
Qee = 4 (0.26 +0.77V4)(.948 €o— Co) (20) 


for (T, -T.) 7 O 


¢ 


Qu 


Qre 


for (ie aS © 


3s ae — V4) (21) 
Du 


) S c er 
on7 Ml ee lg ek 


It is apparent that the contribution from the latter equation 
is small and hence will be considered zero. 

By equating the relationship (19) given by Laevastu for 
sensible heat flux with the integral of equation (18) with respect 
to mass per unit area, we may determine a reasonable value for 
the constant "A". 


Integrating (18) with respect to mass per unit area gives 


ts 


Qe = \AAUIT LLY ede 


0) 


kod 


where 
Introducing the hydrostatic relationship, ie d2= “of and integrating 


with respect to pressure, and then equating this term to the heat 


flux (19) yields 
7, 





Ag Ws\ Bop Ts -Ta) = 34.(02b+ 0.774) (TsTe) Ho) 
ex rr oY ol 4 


Solving this equation for A, and placing £\b)= (0.264, 71V,) 
leads to a value of 63.5 for A. The resulting term for the sensible 


heat contribution becomes 


> ie oe 
Qy = 63.5 (0.ae eet ae) (<8..) (23) 


Furthermore, since the wind uSed in the prognostic model will 
be the geostrophic wind at 1000 mb which is about 130 to 140 per- 
cent greater than the surface wind, the constant 0.77 in equation (23) 
was reduced to 0.5. This constant may be refined to further enhance 
the accuracy of the forecast. 
The final form for this term thus becomes 
Quy = 63.5 oa4 10.84) Ts-Ta)( Ae) (24) 
10 


for ae =e) 7 O 


Bowen (1926) related the convective transfer of sensible heat 


to evaporation by the following relationship 


Qu a SBS = K (Te- Ta.) Pa. (25) 


- (Cs —@a) 1000 


A. H. Gordon al (1952) using meteorological data from 


British ocean weather ships in the North Atlantic related the mean 


18 


values of vapor pressure differences (e,-e.) to individual whole 
degree temperature differences of (Ts-Ta) for wind forces 0-3, 4-6 
and greater than 6. 

Roll (4 using this data computed the Bowen ratio from 
corresponding average temperatures and water vapor differences 
(Ta-Ts), (e,-e,) for wind forces Beaufort 4 and 8. Only a very 
slight increase in the Bowen ratio was observed for force 8 winds. 
The results of those computations are presented in Fig. 2. 

The curve of Fig. 2 was statistically fitted using a linear 
rational transformation and solving for the Bowen ratio in terms of 
temperature difference (Ts-Ta) with the following result: 


= ves = ~A (Ts -To) 
p = | + 21a Calm 2) 


As was previously indicated, the diabatic heating term is to 
be comprised only of the two elements Q, and O. 7 MenCemunc 
use of relationship @6) enables us to eliminate the vapor pressure 
differences as a parameter in the forecast. Thus Q (Total Diabatic 


Heating Rate) = Qs ay Qe 


pe. Sle =1e)Q, (27) 


€ aL (T—Ta) 


iN 


The total heating rate determined from equations (2) and @7) 


is thus 


Q = 311-5 (1 +.34(Ts-TL)) 26 ++SMse)( He) 


(28) 
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From equations (2) and (18), the rate of change of thickness 
due to diabatic heating in the lower and upper layers, which will 


be termed Qo and Qi respectively, may be expressed as follows 
_ RAP (BRYA, - 
= Gee EIAMM)C -T.) 
ee b A 4(\4) (1. - Te) (29) 


QQ) = *I59, Fe AA(Ve )(T Ta) (30) 
or Q, = 5/5 Ca (3) 


We have thus established a simplified distribution of the 


effects of the diabatic heating. 


From equations (28) and (81) the rate of change of thickness 


due to diabatic heating in the lower layer is then 


QQ = ‘ha Q 


= \S.2 (14632 (T~Ta))(abr.5 Visio) (a) 
ize)e (T, eae) 70 


Similarly, it may be shown that 


Q),= Ores (Vasie ) (Te le) (33) 


fon GO 


(32) 


a0 


In comparison with equation 82), it is apparent that the 
magnitude of the thickness change due to this term will be quite 
small. However, the constant .27 may be increased if this 
contribution is deemed too small. 

IV. Development of the terrain pressure and friction term 
According to Berkofski and Bertoni [1] , the vertical velocity 


produced by terrain may be expressed as follows: 


Wy = --4VV £5 (34) 


Making the geostrophic assumption, and also assuming that 
the 1000-mb wind is valid for all levels of terrain, this relationship 


may be expressed as follows: 


W.- = - “ KI (D2, (35) 


It is apparent that this term may be somewhat small in cases 
of high terrain due to the fact that the geostrophic wind at 1000 mb 
is used rather than an integrated wind over the entire height of the 
terrain. However, in nature, the wind will normally tend partially 
to parallel the terrain contours, hence this reduction will somewhat 
compensate for this deficit. 

It was originally intended to use a more refined function for 
terrain induced vertical velocities and to include the effects of 
friction. Time limitations however precluded any further work on 


this term. 


rok 


It is felt that further expansion of this term or the inclusion 
of the terrain term now calculated by FNWF into the program will 
render some improvement to the prognosis. 


V. Programming the Prognostic Equations for the Control 
Data 1604 Computer 


The prognostic equations were put into finite-difference form 
and programmed for the Control Data 1604 computer used in 
conjunction with the FE880 Univac storage drum. The 1604 computer 
is a stored program, general purpose digital computer with a storage 
capacity of 32,768 words in two indpendent storage units. The 
FE880 storage drum is connected to the 1604 computer by a high 
speed channel which conveys information at a rate of 3,000,000 
bits per second. The drum capacity is 256 blocks of 5,120 
48 bit words per block. 

The program was written in symbolic coded relocatable 
assembly program (S CRAP) for a 63 by 63 hemispheric grid, with 
one hourly time steps. 

An option was written into the program to allow for equating 
the vertical velocity.due to terrain and the heating term to zero 
in order to ascertain their respective effects on the prognosis. 

The 750-mb layer is not considered a standard layer in the 
atmosphere, so no synoptic information is available for the 
initialization of the program. To account for this, the following 


equation developed from the thickness equation was derived. 


ae 





om a Og omy 


S775 


— (58.18 (Tew- Too) + 2OAT wo) (36) 


where h is the standard thickness of the layers 700 to 750 mb. 


St7795 
Once the initialization ofthe program is completed the 750 mb 
height may be computed using the forecast thickness for the upper 
level and FNWF's 500-mb barotropic prognosis. 
The absolute vorticities (1, a /4,.) are computed with FNWF's 


subroutine SAR, This routine computes the absolute vorticity over 


a 63 by 63 grid using the following relationships 


h= 4, Sine + Ao (Art Ant As Ay HAs) 


Sie 


Where a, = .0746711 
Ao = .2704348 
To use this routine, all data must be scaled antl The 
absolute vorticity is equated to zero for all points south of 10 N 
latitude. 
After computing the absolute vorticities and reading in the 
various D fields, the Jacobians may be computed utilizing Fleet 


Numerical Weather Facility standard subroutine SAB. The 


Subroutine computes the Jacobian in the following manner 


A = 3 (8,¢) 
= (Ax B)AyC)- (4B) (Ax¢) 


The program inserts zeroes in field A at the boundaries. 


a3 


™ 
The next term computed is the VODs . This term is 
computed by Fleet Numerical Weather Facility's subroutine SAD, 


which computes as follows: 
oy 
A =V"B 
ae Piste 


+L i +B, 7B 


i-\,4 
This program also zeroes the field at the boundaries. 

To compute the heating term, Qo we must obtain an hourly 
surface temperature based on the forecast thickness and also 
compute a 1000-«mb geostrophic wind. The sea surface temperature 
was assumed to remain constant throughout the period of the 
forecast. 

The equation developed to Ae enniné an hourly surface 


temperature as derived from the predicted thickness fields is as 


follows: 


Wee: =34\5 x 10" | (oss h>- 1.33 h) dn (Bans o% Be) +3.4Th| (37) 


After Qo is calculated, each point in the grid is checked to 
determine whether it is a land or sea point. At land points Q5 
is equated to zero; at sea points, as calculated. 

The terrain term Wo is computed uSing the Jacobian subroutine 
and also checked against the land/sea table. In this case, sea 


points are zeroed and land points left unchanged. 


Can 


The terms now calculated are collected and make up the forcing 
function for the Helmholtz equation (14). This equation is solved 
uSing FNWF's subroutine SAE, providing A aa for the first hour of 
the prognosis. 

To solve the Helmholtz aquation, WA-BA =(C, the program scans 
through field A laterally from the lower left point to the upper right 
point using the extrapolated Liebmann method of computing a new 
point value. The boundary conditions have already been determined 
and stored in the A field, the guess field,which is originally set 


to zero. 


V4 y 2, 
VA 


, Be: it — (AB), — Ca a 
14 A. +Al SE ac Ma 1 


¥ = Y 
Roy ai IX Coed _ UNO —C by 4 
* i Cue | 


When for all points [Ri lk € where fis the criterion for 
convergence or the point where the residual R is considered 
liquidated, the field is assumed relaxed and a proper solution has 
been obtained at each grid point. 

When Ah, has been determined, all the parameters required 
for the computation of Ah, are available for use in equation (17). 

The change in height at 1000 mb is then determined by 
subtraction of Ah, plus 4 ho from the 500-mb height change as 


computed by the FNWF barotropic model. 
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6. Results and Conclusions 

Twenty four hour thickness forecasts were produced from 
initial data on 20, 23 and 29 May 1966. The initial program 
contained a number of smoothing passes at many stages of the 
prediction procedure. AS a result the initial thickness prognosis 
was very smooth, almost zonal, as shown in Figure }. Compar- 
ison of this prognosis with the initial and verifying maps, 
Figures 3 and 6, clearly shows that oversmoothing was a major 
source of error. In figures 5, 10 and 13 for the 20th, 23rd 
and 29th of May the results were slightly more encouraging. 
Here smoothing has been reduced somewhat and lower limits 
were placed on absolute vorticity and stability. Still further 
elimination of smoothing throughout the prediction procedure 
gave the thickness and 1000-mb prognosis for the 29th in 
Figures 15 and 19, which show genuine promise as a useable 
forecast. No 1000—-mb analysis existed for the 21st and eth. 

Generally, the pressure systems moved realistically and 
the final positions were good. However, some excess move— 
ment was observed which might be partially Fitribuced to the 
smoothing still remaining in the model. Some overdevelopement 
of both high and low pressure systems in Fiigure 15 was also 
observed. However FE’ NWE' has a filtering subroutine which in 
addition to filtering out the ultra short waves reduces the 
amplitude of these centers. This subroutine diminishes the 


effect of simply spreading the patterns to adjacent grid points 
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which is inherent in smoothing subroutines. Smoothing can also 
give rise to fictitious movement of pressure systems. 

By replacing the smoothing with a filtering pass some of 
the difficulties of excess movement and overdevelopement could 
be overcome. 

Another problem is the tendency of some thermal systems, 
known for their persistance, to move with the general circulation, 
an example is the heat ‘ime in Southwestern United States. 
Periodic comparison with climatology could be fruitful in reducing 
this tendency. 

Due to a lack of available computer time little consideration 
could be given to the diabatic heating and terrain terms. Further 
refinement and inclusion of these terms could result in some 
improvement of the prognosis. 

The program running time on the 1604 computer for a 2h- 
hour forecast was about 30 minutes, including hourly writing 
on magnetic tape of a number of fields which would be unnecessary 
in an operational model. 

The F NWF thickness forecast in Figure 16, which is based 
on simple horizontal advection and numerous constraints, is 
obviously superior to the prognosis resulting from the model 
under investigation... However, the superiority is a result of 
several years of refinement and empirical adjustment. Time 
did not permit such thorough investigation of this research model. 


In conclusion, it has been determined that the assumptions 


a7 


used in the construction of this model are reasonable. With 
proper adjustment and refinement this model should yield a good 
forecast of thickness. This, in turn, should produce a good 


surface prognosis. 


ao 
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Figure 1. Two Layer division of the atmosphere between 
1000—-mb and 500-mb. hy, and hop are the thickness of the 


“layers py to Pr and P19 ol 6 Dy respectively. 
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Bowen Ratio vs. Sea-Air Temperature Difference 


Figure 2. 
(dashed curve) 
Linear rational transformation for positive values, 


Equation (26) (Solid curve). 
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